Stress Echo Descriptive Analysis

Introduction

This report reads an existing R dataset (officially a “data frame”) and makes various descriptive summaries, packaging text, code, and output into an html file. Some of the graphics are rendered both as static and as interactive graphs. The interactive ones allow you to pan, zoom, and hover with the mouse to see more details. Interactive graphics are constructed using the R plotly package, which the Hmisc package loads automatically. To run this report you’ll need to install the Hmisc and plotly packages and their dependencies (RStudio will take care of the dependencies).

Loading User-Contributed R Packages

In R a package is a collection of functions, function documentation (“help files”), and example datasets. When you install R and RStudio many standard packages are automatically installed. There are over 10,000 user-contributed packages that extend what R can do. These include specialty packages for biomedical research (e.g., flow cytometry analysis, genomics, interfaces to REDCap) and a huge number of packages for analysis, table making, and graphics. After using the RStudio Tools menu to install add-on packages, you load packages into R to make them available when running your report by using require or library commands. In the following we get access to the Hmisc package.

Code

require(Hmisc)   # Hmisc must already be installed

If you installed Hmisc after 2022-08-15 you can load this complete report script into your RStudio script editor window by running the command (usually from the R console under RStudio) getRs('stressecho.qmd', put='rstudio').

Accessing Data

There are many ways to import data into R, as described here. Much of the time you will be importing .csv or REDCap files or binary files created by Stata, SPSS, or SAS. For learning R there is a wide variety of ready-to-use datasets on the Vanderbilt Department of Biostatistics dataset repository at hbiostat.org/data. These datasets are already in R binary format and many of them are fully documented and annotated. In the Hmisc package there is a function getHdata that finds these datasets, downloads them from the web server, and loads them into R’s memory for easy access.

Use thegetHdata function to download the stress echocardiogram dataset and read it into R’s memory. This dataset is annotated with variable labels and, for continuous variables, units of measurement. These annotations are used by contents, describe, and certain graphics functions.

Sometimes datasets have long names that are tedious to type. I use a convention of copying the currently active dataset into an R data frame called d to save typing during analysis.

Code

getHdata(stressEcho)
d <- stressEcho   # copy dataset so can refer to short name d

Data Dictionary

The names of variables in the stress echo dataset can be printed by using the command names(d). But let’s use an Hmisc package function to provide more metadata¹

¹ results='asis' in the chunk header prepares R to expect specially formatted html output to be created by the html function.

Code

html(contents(d))

---

Data frame:d

558 observations and 31 variables, maximum # NAs:0  

---

Name	Labels	Units	Levels	Storage
bhr	Basal heart rate	bpm		integer
basebp	Basal blood pressure	mmHg		integer
basedp	Basal Double Product bhr*basebp	bpm*mmHg		integer
pkhr	Peak heart rate	mmHg		integer
sbp	Systolic blood pressure	mmHg		integer
dp	Double product pkhr*sbp	bpm*mmHg		integer
dose	Dose of dobutamine given	mg		integer
maxhr	Maximum heart rate	bpm		integer
pctMphr	Percent maximum predicted heart rate achieved	%		integer
mbp	Maximum blood pressure	mmHg		integer
dpmaxdo	Double product on max dobutamine dose	bpm*mmHg		integer
dobdose	Dobutamine dose at max double product	mg		integer
age	Age	years		integer
gender			2	integer
baseEF	Baseline cardiac ejection fraction	%		integer
dobEF	Ejection fraction on dobutamine	%		integer
chestpain	Chest pain			integer
restwma	Resting wall motion abnormality on echocardiogram			integer
posSE	Positive stress echocardiogram			integer
newMI	New myocardial infarction			integer
newPTCA	Recent angioplasty			integer
newCABG	Recent bypass surgery			integer
death				integer
hxofHT	History of hypertension			integer
hxofDM	History of diabetes			integer
hxofCig	History of smoking		3	integer
hxofMI	History of myocardial infarction			integer
hxofPTCA	History of angioplasty			integer
hxofCABG	History of coronary artery bypass surgery			integer
any.event	Death, newMI, newPTCA, or newCABG			integer
ecg	Baseline electrocardiogram diagnosis		3	integer

---

Category Levels

gender

male

female

hxofCig

heavy

moderate

non-smoker

ecg

normal

equivocal

MI

In the main output you’ll find that the number of levels for categorical variables (R factor variables) are highlighted. These are hyperlinks, and if you click on the number the pointer will jump to the list of levels below.

Descriptive Statistics

The Hmisc package describe function produces appropriate descriptive statistics for each variable in a data frame depending on the variable’s nature (binary, categorical, continuous, etc.). If you just type describe(d) you’ll get plain text output with no graphics. Running the output of describe through the html function converts it to html format which contains small graphics files to render high-resolution spike histograms for continuous variables.

See this for definitions of some of the items in describe output.

Code

des <- describe(d)   # store results of describe() in des
html(des)            # can also use html(describe(d))

d

31 Variables   558 Observations

---

bhr: Basal heart rate bpm

n	missing	distinct	Info	Mean	Gmd	.05	.10	.25	.50	.75	.90	.95
558	0	69	0.999	75.29	16.57	54.0	58.0	64.0	74.0	84.0	95.3	102.0

lowest : 42 44 45 46 47 , highest: 108 115 116 127 210

---

basebp: Basal blood pressure mmHg

n	missing	distinct	Info	Mean	Gmd	.05	.10	.25	.50	.75	.90	.95
558	0	94	0.998	135.3	23.35	104.0	110.0	120.0	133.0	150.0	162.3	170.1

lowest : 85 88 90 97 98 , highest: 192 194 195 201 203

---

basedp: Basal Double Product bhr*basebp bpm*mmHg

n	missing	distinct	Info	Mean	Gmd	.05	.10	.25	.50	.75	.90	.95
558	0	441	1	10181	2813	6607	7200	8400	9792	11663	13610	14770

lowest : 5000 5220 5280 5400 5460 , highest: 17604 17710 17748 21082 27300

---

pkhr: Peak heart rate mmHg

n	missing	distinct	Info	Mean	Gmd	.05	.10	.25	.50	.75	.90	.95
558	0	105	1	120.6	25.36	81.85	90.70	106.25	122.00	135.00	147.00	155.15

lowest : 52 61 62 63 66 , highest: 170 171 176 182 210

---

sbp: Systolic blood pressure mmHg

n	missing	distinct	Info	Mean	Gmd	.05	.10	.25	.50	.75	.90	.95
558	0	142	1	146.9	40.72	96	102	120	141	170	200	210

lowest : 40 60 70 79 80 , highest: 240 250 274 283 309

---

dp: Double product pkhr*sbp bpm*mmHg

n	missing	distinct	Info	Mean	Gmd	.05	.10	.25	.50	.75	.90	.95
558	0	508	1	17634	5765	10256	11341	14033	17060	20644	24536	26637

lowest : 5100 5940 7490 8100 8360 , highest: 32518 33400 33840 38205 45114

---

dose: Dose of dobutamine given mg

n	missing	distinct	Info	Mean	Gmd
558	0	7	0.84	33.75	8.334

lowest : 10 15 20 25 30 , highest: 20 25 30 35 40

 Value         10    15    20    25    30    35    40
 Frequency      2    28    47    56    64    61   300
 Proportion 0.004 0.050 0.084 0.100 0.115 0.109 0.538
 

---

maxhr: Maximum heart rate bpm

n	missing	distinct	Info	Mean	Gmd	.05	.10	.25	.50	.75	.90	.95
558	0	103	1	119.4	24.64	82.0	91.0	104.2	120.0	133.0	146.0	154.1

lowest : 58 62 63 66 67 , highest: 170 171 176 182 200

---

pctMphr: Percent maximum predicted heart rate achieved %

n	missing	distinct	Info	Mean	Gmd	.05	.10	.25	.50	.75	.90	.95
558	0	78	0.999	78.57	16.86	53	60	69	78	88	97	104

lowest : 38 39 40 41 42 , highest: 116 117 126 132 133

---

mbp: Maximum blood pressure mmHg

n	missing	distinct	Info	Mean	Gmd	.05	.10	.25	.50	.75	.90	.95
558	0	132	0.999	156	35.03	110.0	120.0	133.2	150.0	175.8	200.0	211.1

lowest : 84 90 92 93 96 , highest: 240 250 274 283 309

---

dpmaxdo: Double product on max dobutamine dose bpm*mmHg

n	missing	distinct	Info	Mean	Gmd	.05	.10	.25	.50	.75	.90	.95
558	0	484	1	18550	5385	11346	12865	15260	18118	21239	24893	27477

lowest : 7130 8100 8360 9240 9280 , highest: 32518 33400 33840 38205 45114

---

dobdose: Dobutamine dose at max double product mg

n	missing	distinct	Info	Mean	Gmd
558	0	8	0.941	30.24	10.55

lowest : 5 10 15 20 25 , highest: 20 25 30 35 40

 Value          5    10    15    20    25    30    35    40
 Frequency      7     7    55    73    71    78    62   205
 Proportion 0.013 0.013 0.099 0.131 0.127 0.140 0.111 0.367
 

---

age: Age years

n	missing	distinct	Info	Mean	Gmd	.05	.10	.25	.50	.75	.90	.95
558	0	62	0.999	67.34	13.41	46.85	51.00	60.00	69.00	75.00	82.00	85.00

lowest : 26 28 29 30 33 , highest: 89 90 91 92 93

---

gender

n	missing	distinct
558	0	2

 Value        male female
 Frequency     220    338
 Proportion  0.394  0.606
 

---

baseEF: Baseline cardiac ejection fraction %

n	missing	distinct	Info	Mean	Gmd	.05	.10	.25	.50	.75	.90	.95
558	0	54	0.994	55.6	10.71	32	40	52	57	62	65	66

lowest : 20 21 22 23 25 , highest: 74 75 77 79 83

---

dobEF: Ejection fraction on dobutamine %

n	missing	distinct	Info	Mean	Gmd	.05	.10	.25	.50	.75	.90	.95
558	0	60	0.992	65.24	12.38	40.0	49.7	62.0	67.0	73.0	76.0	80.0

lowest : 23 25 26 27 28 , highest: 86 87 89 90 94

---

chestpain: Chest pain

n	missing	distinct	Info	Sum	Mean	Gmd
558	0	2	0.64	172	0.3082	0.4272

---

restwma: Resting wall motion abnormality on echocardiogram

n	missing	distinct	Info	Sum	Mean	Gmd
558	0	2	0.745	257	0.4606	0.4978

---

posSE: Positive stress echocardiogram

n	missing	distinct	Info	Sum	Mean	Gmd
558	0	2	0.553	136	0.2437	0.3693

---

newMI: New myocardial infarction

n	missing	distinct	Info	Sum	Mean	Gmd
558	0	2	0.143	28	0.05018	0.09549

---

newPTCA: Recent angioplasty

n	missing	distinct	Info	Sum	Mean	Gmd
558	0	2	0.138	27	0.04839	0.09226

---

newCABG: Recent bypass surgery

n	missing	distinct	Info	Sum	Mean	Gmd
558	0	2	0.167	33	0.05914	0.1115

---

death

n	missing	distinct	Info	Sum	Mean	Gmd
558	0	2	0.123	24	0.04301	0.08247

---

hxofHT: History of hypertension

n	missing	distinct	Info	Sum	Mean	Gmd
558	0	2	0.625	393	0.7043	0.4173

---

hxofDM: History of diabetes

n	missing	distinct	Info	Sum	Mean	Gmd
558	0	2	0.699	206	0.3692	0.4666

---

hxofCig: History of smoking

n	missing	distinct
558	0	3

 Value           heavy   moderate non-smoker
 Frequency         122        138        298
 Proportion      0.219      0.247      0.534
 

---

hxofMI: History of myocardial infarction

n	missing	distinct	Info	Sum	Mean	Gmd
558	0	2	0.599	154	0.276	0.4004

---

hxofPTCA: History of angioplasty

n	missing	distinct	Info	Sum	Mean	Gmd
558	0	2	0.204	41	0.07348	0.1364

---

hxofCABG: History of coronary artery bypass surgery

n	missing	distinct	Info	Sum	Mean	Gmd
558	0	2	0.399	88	0.1577	0.2661

---

any.event: Death, newMI, newPTCA, or newCABG

n	missing	distinct	Info	Sum	Mean	Gmd
558	0	2	0.402	89	0.1595	0.2686

---

ecg: Baseline electrocardiogram diagnosis

n	missing	distinct
558	0	3

 Value         normal equivocal        MI
 Frequency        311       176        71
 Proportion     0.557     0.315     0.127
 

---

describe also has a plot method that translates descriptive statistics to purely graphical form and draws larger high-resolution histograms for continuous variables. The plot method produces two graphics objects: one for categorical variables and one for continuous variables. We save the overall plot result in object p then reference two sub-objects²

² For these static graphics we could have just typed plot(des) and two plots would appear in succession. That approach does not work when making interactive graphs later.

Code

p <- plot(des)

Code

p$Categorical

Figure 1: Categorical variables

Code

p$Continuous

Figure 2: Continuous variables

To render these two graphs in interactive format using plotly we set a system option that the Hmisc package monitors: grType. It’s default value is 'plain'.

Code

options(grType='plotly')
p <- plot(des)

Code

p$Categorical

Figure 3: Categorical variables. Hover over points to see numerators and denominators of proportions, and hover over the leftmost point to see more information.

Code

p$Continuous

Figure 4: Continuous variables. Hover over the leftmost spike to see a full list of descriptive statistics. Hover over other spikes to see \(x\)-axis values and frequency counts.

Computing Environment

A function session stored in an object in the Hmisc package prints the current computing environment in a nice format.

 R version 4.2.1 (2022-06-23)
 Platform: x86_64-pc-linux-gnu (64-bit)
 Running under: Pop!_OS 22.04 LTS
 
 Matrix products: default
 BLAS:   /usr/lib/x86_64-linux-gnu/blas/libblas.so.3.10.0
 LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.10.0
 
 attached base packages:
 [1] stats     graphics  grDevices utils     datasets  methods   base     
 
 other attached packages:
 [1] Hmisc_4.7-2     ggplot2_3.3.5   Formula_1.2-4   survival_3.4-0 
 [5] lattice_0.20-45
 

To cite R in publications use:

R Core Team (2022). R: A Language and Environment for Statistical Computing. R Foundation for Statistical Computing, Vienna, Austria. https://www.R-project.org/.

To cite the Hmisc package in publications use:

Harrell Jr F (2022). Hmisc: Harrell Miscellaneous. R package version 4.7-2, https://hbiostat.org/R/Hmisc/.

To cite the ggplot2 package in publications use:

Wickham H (2016). ggplot2: Elegant Graphics for Data Analysis. Springer-Verlag New York. ISBN 978-3-319-24277-4, https://ggplot2.tidyverse.org.

To cite the survival package in publications use:

Therneau T (2022). A Package for Survival Analysis in R. R package version 3.4-0, https://CRAN.R-project.org/package=survival.